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We set up a harmonic lattice model for 2D defect melting which, in contrast to earlier simple-cubic 
models, lives on a triangular lattice. Integer- valued plastic defect gauge fields allow for the thermal 
generation of dislocations and disclinations. The model produces universal formulas for the melting 
temperature expressed in terms of the elastic constants, which are different from those derived 
for square lattices. They determine a Lindemann-like parameter for two-dimensional melting. In 
contrast to the square crystal which underwent a first-order melting transition, the triangular model 
melts in two steps. Our results are applied to the melting of Lennard- Jones and electron lattices. 
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I. INTRODUCTION 

Melting transitions Q, |^ |^ are of both technical and 
theoretical interest. Under suitable experimental circum- 
stances one can study the melting process in two dimen- 

[ sions (2D). If the interparticlc interactions are simple, the 
lowest-energy crystalline order in 2D is a close-packed tri- 

' angular. It is observable for electrons on a liquid-helium 
surface for adsorbed noble-gas atoms floating on an 
incommensurate substrate 0, for colloidal suspensions 
of highly charged spheres confined between parallel glass 
plates for paramagnetic colloidal crystals JJ, and for 

] vortex lattices in thin films .8J. Simulation experiments 
provide us with further examples: rigid disks Gaus- 

I sian cores |lOj . and particles interacting with Lennard- 
Jones as well as Coulomb forces (Wigner lattice) ^3 ■ 

I There are various theories of 2D melting 0, |3i Q> Qi 
[T5I [T^ . ITsj . Most popular is the phenomenological 
theory of Halperin, Nelson, and Young 14] which, in- 
spired by the Kosterlitz and Thouless (KT) theory of 
vortices in superfluid films explains the transition 
by the statistical behavior of defects. This KTHNY the- 
ory suggests that melting could proceed in a sequence 

' of two continuous KT transitions, the first caused by the 
unbinding of dislocations, the second by the unbinding of 
disclinations. Such a sequence was indeed found in some 
2D crystals experimentally and by simulations. 

There are also explicit models on square lattices which 
allow to study defect melting in detail. These models are 
minimal in the sense that they contain the correct har- 
monic elastic energy and no anharmonic terms. The de- 
fects arise from integer- valued plastic gauge fields which 
account for the thermal creation and annihilation of dislo- 
cations and disclinations. The simplest such model shows 
only a single first-order melting transition |l7j| . Only 
after introducing an extra higher-gradient elastic term 
which gives the crystal a tunable angular stiffness, the 
first-order transition separates into two successive KT 
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transitions [iSj. This splitting can be observed in the 
laboratory, for instance in Xenon overlayer lattices on 
graphite pOl |. 

Presently investigated systems live mostly on trian- 
gular lattices, and these seem to melt in two steps 
Ellll^lSlEllIlElll^, except for the Gaussian core 
model. The question arises whether a minimal melting 
model living on a triangular lattice will always undergo 
two successive transitions at the level of first gradient 
elastic energies only, so that no extra higher-gradient 
terms are needed for splitting them as on square lattices. 
It is the purpose of this paper to answer this question to 
the affirmative. We shall generalize the previously con- 
structed lattice models from square lattices H, [ill 
to triangular lattices and find, indeed, strong indications 
that they melt in two KT transitions. 

A universal feature of all such defect models is that it 
is the combined separation of dislocations and disclina- 
tions which makes the melting transition first-order |l7j . 
Dislocations alone would cause a similar phase transition 
as in superfiuid helium, which undergoes a continuous 
transition. 

An important virtue of the lattice defect models on 
square lattices was that they lead to a simple univer- 
sal melting formula determining the melting point in 
terms of the elastic constants. The result is obtained 
from a lowest-order approximation, which combines the 
high-temperature expansion of the defect contributions 
to the free energy density with the low-temperature ex- 
pansion. The melting temperature is determined by the 
intersection of the two curves. This turns out to produce 
precisely the time-honored Lindemann criterion. In addi- 
tion, however, there is a prediction of the universal value 
of the Lindemann number where melting occurs. The ac- 
curacy of the approximation was demonstrated quantita- 
tively for square lattices in 3D and 2D 0, . Recently, 
the results were successfully extended to face-centered 
and body-centered cubic lattices in 3D ^2^ . We shall find 
a similar formula for triangular lattices in 2D. The mini- 
mal lattice defect model will be constructed in Section II. 
From its partition function we calculate in Section III the 
melting temperature. In Section IV, we compare our re- 
sult with the melting transition of Lennard-Jones and 
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electron lattices. Section V discusses the melting values 
when taking into account nonlinear elasticity effects. In 
Section VI we use the model to predict a universal gen- 
eralized Lindemann number for 2D lattices. 



II. LATTICE HAMILTONIAN WITH DEFECTS 

The elastic energy of a crystal in the continuum ap- 
proximation is given in D dimensions by (28| 



(1) 



where repeated indices are summed. The elastic constant 
H is the shear modulus, and A the Lame constant. In 
two dimensions, these constants are related to the elastic 
constants Cy with lattice symmetry by /i = Cge and A = 
cii — 2c66. The combination cn — cqq is the modulus of 
compression. The fields Ui{x) {i = x,y) are the cartesian 
components of the atomic displacements. 

In a first step we construct a Hamiltonian on a 2D 
triangular lattice whose continuum elastic energy agrees 
with (Q. Of course, this construction is not unique. We 
restrict the freedom by the following extension require- 
ments: 

1. The lattice derivatives V^; and Vj, must have the 
continuum limits dx and dy, respectively. 

2. V^^, Vj, should be maximally symmetric linear com- 
binations of the nearest neighbor lattice derivatives 
V(;), / = 1,2,3 with respect to the point group of 
the lattice. 



3. V^: and Vj, should have the same transformation 
properties as dx and dy under the action of the 
point group of the lattice. 



Eq. Q, but with the replacement e(j) — s- ej^). 
definition follows the identity 

(V(i)+V(2)+V(3))/(x)==0 



From the 



(4) 



for any function /(x) on the lattice sites. 

The 2D point group around each lattice face is given 
by Csv 11^' This group causes permutations of the link 
directions e(;) and the derivates V(;). These can be rep- 
resented by matrices which form a 3 x 3-dimensional 
representation of the point group C^^ in the space of 
linear combinations of the link vectors or of the lattice 
derivates. 

We now define cartesian lattice versions V^; and Vj, 
of the derivatives dx and dy as C^y- symmetric linear 
combinations of the lattice derivatives V(;): 



- e{l)x,y^{l) 



V x,y — 2 

Inserting the components from these become 



(5) 



dx 



Vx = 




1 1 

-V(i)-V(2) + -V(3) 



V3. 



V3. 



— V(3) 



(6) 



It is easy to check that these fulfill the above require- 
ments 1-3. 

Without proof we mention that by using group theo- 
retical methods (see e.g. we can show that the ex- 
tension requirements 1 — 3 are sufficient to obtain unique 
lattice derivates V^; and Vy. This is true for all other lat- 
tice symmetries. The only property needed for the proof 
is that the lattice possesses a link vector which generates 
all others by applying the symmetry group of the lattice. 



III. LATTICE HAMILTONIAN 



In the following we denote by e^;) the oriented link vec- 
tors of the triangular lattice which surround a triangular 
face of the lattice. There are three such vectors: 



6(1) 
6(3) 



277 . 277 

cos — , sm — 
6 6 



3(2) 



-1,0), 



2n 



2tt 



cos — , — sm — 
6 ' 6 



(2) 



and three conjugate vectors e(;) = — e(;) 

For each of the link vectors there exist a lattice deriva- 
tive defined by 

V(i)/(x) = [/(x + e(i))-/(x)]/a, 
V(2)/(x) = [/(x)-/(x-e(2))]/a, 
V(3)/(x) = [/(x-e(2))-/(x + e(i))]/a. (3) 

where a is the lattice spacing. In addition, there is a 
conjugate derivative V(i)/(x) defined by — V(i)/(x) as in 



It is now easy to set up the lattice version of the con- 
tinuum Hamiltonian (Q. We simply replace the partial 
derivatives by the lattice derivatives lO, and obtain 

X 

(7) 

where the sum runs over all lattice sites, and v = 
\/3a^/2 is the area of the fundamental cell of the lattice 
(for the square lattice v = a^). As in the earlier work on 
square lattices 2] we shall first consider a symmetrized 
version of the Hamiltonian {Tj) which has the advantage 
of leading to explicit formulas for the melting transition 
for all /i, A. This arises by replacing by in the A 
term: 

Ei^t - {|[V.u,(x) + V,u,(a;)]2 + A[V,u,(x)]2} . 

X 

(8) 
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Without an explicite discussion in this paper we mention 
that one can show the difference in the results using the 
symmetrized Hamiltonian JHJ and the Hamiltonian 
to be negligible for triangular as well as square lattices. 
The most symmetric way of writing ((7J is 

mat— — -. — 

4 ' — ' I V] \ ' — ' z 



H j ^ ( H -\Y1 "^{^Hi)] 



\ kij / 



(9) 



Here we have expanded the displacement vector in (j?)) in 
a symmetric way to u(a;) — U(^i^{'x)e(j_y This representa- 
tion of u(x) is not unique due to the overcompleteness of 
the three link vectors It becomes unique by setting 
M(3) = and taking into account yielding 

Eut = ^- 2^\Cll[ V(1)U(1) + V(2)U(2) 



+ g C66 |^V(i)U(i) - V(2)W(2) + 2V(2)U(i) - 2V(i)-U(2) j 

+ 4C66^V(1)U(2)V(2)M(1) - V(i)M(i)V(2)U(2)^ > • (10) 



The classical partition functions of the triangular crystal 
including dislocation and disclination degrees of freedom 
is then given by 



Z 



n 



dMj(x) 



E 

ii™(x)} 



'^[n(im)\ exp[-£;iat/fcBT]. 

For simplicity, we choose periodic boundary conditions 
for the displacement fields u(x). Since defects are 
spanned over an integer number of fundamental cells 
of the lattice, the functional $[rt(;TO)] should restrict 
the range of the integer- valued gauge field Fur- 
thermore, it acts also as a gauge-fixing, in such a way 
that only physically independent degrees of freedom are 
summed. 

Due to the identity we have the freedom to carry 
out the replacement 'S/x,y — > '^x,y + a2;,j/(V(i) -I- V(2) -t- 
V(3)) in © with some constants ax,y without changing 



the elastic energy in JTj) or (|8|) . One can see that by con- 
sidering relation (|ll|l , the resulting Hamiltonian and also 
the partition functions depends strongly on ax.y This is 
because the replacement (|ll|l does not take into account 
the fact that defect lines as dislocations and disclinations 
are built of interstitials or vacancies, respectively, which 
allows only one interstitial or vacancy per fundamental 
cell. This restricts the Volterra surfaces to run through 
an integer number of fundamental cells. (See Fig. 1). 
This restriction on the defect gauge fields results in the 
following constraint 



(x) + n(2/) (x) + n(3,) (x) = 0, 



(13) 



IV. INCLUDING DEFECT GAUGE FIELDS 

With the goal of studying defect-induced melting tran- 
sitions, the lattice representations ((TJ and ((Sjl of lin- 
ear elasticity must be extended by integer-valued defect 
gauge fields, the lattice version of plastic fields. From 
the textbook 0] we know how to do this for square lat- 
tices. The defect gauge fields reflect the fact that due to 
fluctuations, atoms are capable of exchanging positions 
with their neighbors and migrate eventually through the 
entire crystal. This process of self-diffusion makes it im- 
possible to specify the displacement fleld uniquely. Thus, 
as a matter of principle, the displacement fleld is multi- 
valued. It is determined only up to an arbitrary lattice 
vector, since it is impossible to say whether an atom is 
displaced by u(x) or by u(x)-f aAf(/)(x)e(/) where iV(;)(x) 
is integer. The displacement vectors may have jumps 
across so-called Volterra surfaces which are described by 
integer- valued defect gauge flelds n(;m)(x) in the lattice 
model 2]. These are, in general, not symmetric. 

The defect gauge fields are included into the lattice 
derivatives of the Hamiltonian ^ by the replacement 

2 



Va;,j^u(x) 



for all / — 1,2,3. Here we take into account that for 
a dislocation the jump over the Volterra surface is 
constant. Disclinations are built of dislocations within 
our lattice model. By taking into account this constraint 
on the defect gauge fields the resulting partition function 
does not longer depend on the replacement constants 



We now must identify and eliminate the gauge degrees 
of freedom in Eq. ()12|l. The overcompleteness of the 
three basis vectors ej-;) of the two-dimensional lattice im- 
plies that the replacement leads to the same de- 
fects if the integer numbers satisfy n(;i)(x)e(i) — for 
I — 1,2,3. The rest of the gauge freedoms can be elimi- 
nated by analogy with the gauge fixing on the square lat- 
tice in the textbook Initially, the displacement fields 
cover only a single fundamental cell. We can extend the 
range to the entire crystal by fixing two of the last four 
components of the defect gauge field rui^y One more 
component can be eliminated since i?iat depends on the 
derivates "S/ xUy and VyUx via the sum "S/ yUx + "S/ xUy- 

Finally, we must take care of the constraints due to 
the periodic boundary conditions. The technical aspects 
of the elimination of the gauge degrees of freedom are 
deferred to Appendix A. Eventually we arrive at the par- 
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FIG. 1: Cut in a triangular lattice. For each lattice site, the 
shaded regions show the faces which are surrounded by the 
oriented link vectors ||5J associated with the lattice derivates 
V(;) in Eq. 0. Dashed line indicates a Volterra surface cross- 
ing an integer number of fundamental cells. 



tition function: 



Wi(x) 



E 



i(x) = 




X I (1 - <5x,b) + <5x,H 

"(ll)(x):n(2l)(x) = -oo 



with the lattice energy ((TJ and the replacements of the 
lattice derivatives 



(V:,,Vy) X u(x) ^ -(e(,)^,e(,)^) x V(,)u(x) (15) 

; 

2 V +73^(21) (x) -2n(ii)(x) - n(2i)(x) 

where B contains the sites at the rightmost and top 
boundary of the sample. 

The gauge degrees of freedoms are fixed by choosing 
n(x) = for X S S, further Va;n(x) = on the rightmost 
boundary Br, and Vyn(x) = on the topmost bound- 
ary Bt of the model. This will be summarized by the 

statement 



n(x) =0 for xe B' . 



(16) 



where B' is the extended boundary consisting on the 
sites B, Bf + ae(i) and Br + ae(2). The boundary fields 
n(]^x') (x)and n(-2i)(x) in (|15|) result from the periodic 
boundary conditions of the displacement field u, which 
restrict the gauge degrees of freedom. 

The final partition function l|14|) describes a triangular 
crystal with harmonic clastic fluctuations and fluctuating 
dislocations and disclinations. 

By a standard duality transformation it is possible to 
rewrite this partition function in a canonical form. Intro- 
ducing the conjugate stress fields crij(x) via an auxiliary 



integral, we see that (|14|l is equal to 

3N/2 



z = 





N/2 ^ 


[4(A + m)J 





A 



n 

X,l<J 



n 



00 



E 

n(x) — — 00 



rfu(x) 



exp ■ 



-y 



2(3 



E4 + JE^' 



A 



4(A + m) 



(E-)^ 



X exp 



[i)Ui ai; 



(17) 



(e(i)^V(;)U2 + e(,)2V(;)Ui) (712 



where /3a = aVV3/2fcBT(27r)2. For details see 0. 

The energy H^n couples the stress field (7^ to the in- 
teger value fields n(y)(x): 



H„n = i 27r ^ — CT12 n(x) + In 



n 



E 



xeB Ln(ii),n(2i)(x)=-oo. 



exp 



j27rQcr(n)(x)n(2i)(x) - ^(7(22) (x) (n(2i)(x)-|- 



-2n(ii)(x)) + -j^ cr(i2)(x) (n(2i)(x) - n(n)(x)) 



(18) 



This is the dual representation of (|14(l . For very large lat- 
tices, the integer- valued fields (x), n(2i) (x) have, of 
course, no physical implications since they exists only on 
the boundary of the system. For the present calculations 
on a finite lattice they are, however, needed. 

The partition function H17|) corresponds to the Hamil- 
tonian ((TJ including defects. The partition function for 
the symmetrized Hamiltonian (jSJ is given by the same 
expression with the replacement 



A 



4(A + /.) 



E^' 



A 



4(A + m) 



(19) 



in the third term of the first exponent in H17|) . In the 
following sections, we shall use the symmetrized Hamil- 
tonian corresponding to the partition function 1171) 
with 1(1?^ as the basis for calculating physical quantities. 



V. MELTING TRANSITION 

We now determine the melting temperature of the 
model. 
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A. Lowest-Order Approximation 



dent components are integrated out using the formulas 



It has been shown in the textbook ^] that the inter- 
section of the free energies of the two hmiting curves of 
the high- and low- temperature expansion yields good 
estimates for the melting temperatures for simple cubic 
lattices in three dimensions and for square lattices in two 
dimensions (see Figs. 12.1 and 12.2 on pp. 1084-1085). 

For small T, the defects are frozen out, and the parti- 
tion function has the classical limit of the partition func- 
tion of the Hamiltonian © with the derivate substitution 
(|TT]l and n(x) = 0: 

Zr-..^(2,,.)-'[j^) e'""-^ (20) 

The parameter i/\ denotes the trace of the logarithm of 
the triangular lattice Laplacian divided by the number 
of sites. The Fourier transforms of the lattice derivates 
© are 



1 



^(1) = -(e^"^""'^' - 1) , 
^ ' ia 



(2) = —(1 - e 

^ ' la 



-iak-e 



•6(2) ^ 



(21) 



The conjugate lattice derivatives have the complex- 
conjugate Fourier transforms A'(/) = ^(i)- Using these, 
we calculate from the momentum integral 



1 



2Abz jbz 
1.22. 



d^fclog 



(22) 



More explicitly, the argument of the logarithm is 



(23) 



4 (cos ak • efx") + cos ak • e^2) + cos ak ■ ©(3)) 

3 



The momentum integral in H23fl runs over the 2D 
BriouUin zone whose area is Abz = (27r)2v^/2a2. Let 
us compare the value J22() with the corresponding on the 
square lattice in Ref. 2] , where the low-temperature par- 
tition function looks like H20I) . but with £a ~ 1-22 re- 
placed by In ~ 1.14, and /3a = ^i^^i/2kBT{2^:Y by 

Above the melting point, the partition function is cal- 
culated from the dual representation (|17() . In the high- 
temperature limit, the defects are prolific. The sum over 
n(x) enforces cri2(x) to have discrete values such that 
terms with cri2(x) ^ are exponentially small for large 
temperature. Thus we can restrict ourselves to terms 
with (712 (x) = only. After carrying out the integrals 
over the lattice displacements u(x), we obtain the con- 
straint that the stress field is divergenceless. Its indepen- 



n 



= exp 
= exp 
and 



dan 



N 



^ [ ye(i)^V(,)(7,i(x) 



(24) 



2Abz Jbz 
N 



2A 



/ d^fc log (2 — 2 cos ak • 6(2)) 
BZ Jbz 



n 



da' 



22 



^ [ ye(,)^V(i)a,;2(x) 



(25) 



= exp 
= exp 



N 



d k log 



- (6 — 4cosak • ei 
3 ^ 



2^BZ Jbz 
— 4 cos ak • 63 + 2 cos ak • 62) 

The partition function Z in this limit is 

N/2 



1.15 



N 



Zt^oo = (27r/3A)"'^/' 
where Ca is the constant 



4(A + p) 



Ca w 0.57. 



CT , (26) 



(27) 



From the intersection of the high-temperature expan- 
sion (|26|l with the the low-temperature expansion (|20|l we 
obtain the lowest-order result for the melting transition 
of a triangular lattice 



Pi^{\ + v) « — e^^^CA «0.52 
47r 



(28) 



where (1 -I- i^) = 2(A + -f- 2[i). This is somewhat 

smaller than /3n(l + j^) ~ 0.81 for the square lattice .2]. 

Let us now see how these results are changed by a sys- 
tematic improvement of the lowest-order low- and high- 
temperature expansions H20|l and (|26|l . 



B. Beyond the Lowest Order 

In the low-temperature regime, we include the lead- 
ing contributions with nonzero defect gauge fields ri(x) 
in the partition function (|14|l . In the dual representa- 
tion (|17() . we must take into account the integer nature 
of n(x). This was done in Ref. Q by reexpressing H17|l 
in terms of the physical defect density field, which is the 
double-curl of the defect gauge field. Here we carry out 
a similar calculation for the triangular lattice. After per- 
forming, in 117|l . the integrals over the displacement field 
u(x) and the stress fields oij (x) we obtain the improved 
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low-temperature approximation Z = Zt^qZ^^ with the 
defect correction factor 



X 7l(x) — — OO 



X exp 



- 47r^/3A(l + i^) E V(i)e(i)a;V(„)e(™)j; n(x) 



"(l2),n(ll)(x) = -ooy 



x,x' 



4a' 



X i;a(x-x') — V(fe)e(fe)^V(„)e(„)yn(x') 



(29) 




The numbers n(x) fulfill the boundary condition In 
this way, we have expressed the low-temperature correc- 
tions to the harmonic partition function as a partition 
function of integer-valued fields n(x) in the bulk and two 
integer valued fields n(ii)(x), n(2i)(x) on the boundary 
of the sample. 

The interaction potential ua is given by the inverse 
square of the triangular Laplacian 



wa(x) 



4a^ 



-VmV 



(/) V(„i)e(/) • e(m) 



(32) 

It is not easy to obtain the leading-order contributions 
of Z^^. This is mainly due to cancellation effects in 
the configuration sum coming from the interaction of the 
fields ?T.(x) and n(]^x')(x), ri(2i)(x) on the boundary of the 
sample. A more efhcient calculation is possible if we go 
over from the defect gauge fields to the defect density 
field given by 

Vi^) = — efc;emn VfcV,„ e(i)^ "•(ii)(x) e^-)^^ . (33) 

The defect field 77(x) describes disclination degrees of 
freedoms 0. The dislocations arise from dipoles of two 
nearby disclinations. The main advantage of ?y(x) is that 
it is invariant under defect gauge transformations of the 
integer-valued fields Thus, we obtain for Z"^^ of 

Eq. 123): 

= ]J ^ g-47r"/3A(l + '>)»)(x)DA(x-x')j)(x) ^ ^g^^ 
X ,,(x) 

The sum X]77(x) ^'uns over all different field configurations 
77 (x) defined by H33|) . The values for wa(x) can be calcu- 
lated via the Fourier transform 

exp(ikx) 



wa(x) 



d'k- 



(35) 




FIG. 2: Basic configurations of defect density r;(x) associated 
with the defect gauge field n(ij)(x) = (5i,i(5j,2 5x,xo (a), and 
n(i_,)(x) — 5i,i5j,i5x,xo (b). The numbers in the triangles 
denote the multiplicity of basic charges 2/3. 




FIG. 3: Defect configuration ?7(x) for n(ij)(x) = ((5i, 2(5^,1 — 
5i,i5j,2 — Si,i5j^i) 5x,xo- Its contribution to the free energy is 
given by the first term in Eq. 1361 . 



We fix the gauge freedom due to overcounting by set- 
ting n;3 — 0. Upon further using H13|l we can ex- 
press 77(x) in terms of the fields ri(jj) with i,j = 1,2 
only. The explicite relations are given in Eqs. (|A1|) and 
(|A2|) . Figure 2(a) illustrates the non-zero 77 (x) values 
for n(ij)(x) — i5i,i(5j,2 ^x.xo: while 2(b) does the same 
thing for n(jj)(x) = (5x,xo- The cases n(y')(x) = 

(5^2 <5x,xo and n(ij)(x) = Si^2 (5^,2 ^x,xo follow from 
linear combinations of 2(a) and 2(b) after a clockwise 
rotation by an angle 27r/3. The space of all 77-fields 
are spanned by linear combinations of these four con- 
figurations modulo translations. The figure shows that 
the values of the ?7-field are multiples of 2/3 which we 
call defect charges. All ?7-configurations are composed 
of two neutral rhombuses which are illustrated by the 
dashed and straight lines in Fig 2. Besides being neutral, 
the rhombuses have also no dipole moments. The fig- 
ure shows the multiplicity of the basic charges for both 
rhombuses. The fact that the //-fields are neutral and 
dipole-free agrees with the situation on the square lat- 
tice. There is one difference, however. On the triangular 
lattice there are neutral dipole-free charge configurations 
composed of multiples of basic charges 2/3 which can- 
not be obtained from an T^-field. One simple example 
is a single basic rhombus shown by the dashed line in 
Fig. 2(a) built of defect charges ±2/3. The complica- 
tion with respect to the square lattice where the Ty-fields 
are built of all neutral and dipole-free charge configura- 
tions is due to the constraint (|13() . Finally we mention 
that the basic configuration h{x) — (2/V3) Sx.xo in H3U|) 
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with ?7(x) 



^Va;Vj,n(x) for x ^ S corresponds to the n(ij)(x). The result is 



defect configuration 2(a). 

We shall now prove that all localized ?7(x) fields can 
be built from localized 'T-(ij) (x)-fields if the sample is in- 
finite. This can be seen from the fact that the line B in 
(I29|l and (|30|l can be moved also to the middle of the sys- 
tem where we held fixed the boundary B' in (|16() . This 
follows immediately from the considerations in Appendix 
A and the periodic boundary condition of the displace- 
ment fields Mi(x) used here. Thus, we have to take into 
account only a linear combination of the most localized 
fields ) (x) in H33() to get the leading order in the defect 
partition function H34|l . 

Writing Z/x as an exponential of the connected dia- 
grams Z = e~^^ , we find the free energy of the defect 
fields 

... (3^) 

where V^, Vy are the triangular lattice derivates in © 
and 



VA (x - 



(37) 



By carrying out a numerical integration of the Fourier 
transform similar to (|22|l we obtain 



{)a(0) w 0.163 , i)A(ae(i)) 



-0.03. 



(38) 



These numbers are found by numerical integration of the 
Fourier representation similar to 1)35(1 . The first term in 
H36|l corresponds to the ?7(x)-configuration in Fig. 3. The 
second term in (|36|l corresponds to the ?7-configuration 
2(a) in Fig. 2. The factor 6iV comes from the six possi- 
bilities to cover this basic 7/-configuration on the lattice, 
this means the configuration 2(a), the rotated configu- 
ration by angle 27r/3 and An/'i plus the negative of all 
these three configurations. The other basic diagrams are 
approximately a factor exp[— 47r^WA (0)/3a (1 + v)] smaller 
than the leading terms in (|36|) . When using the mean- 
field value /3a (1 + v) ~ 0.6 we obtain additive correc- 
tions which approximatively a factor 0.04 smaller than 
the leading terms in — F^^^/cbT calculated above. From 
this, we derive also a temperature regime for which our 
lowest defect configuration result can be trusted. This is 
given by the inequality 47r^i;A(0)/3A(l + v) > I resulting 
in Pa{1 + v) > 0.15. 

Recall that for square lattices, the defect correction 
factor looks similar to H36|l. except that the second term 
is missing, and that wa in the first term is replaced by 
where iP is defined by Eq. (|37|l . with being the 
inverse square of the Euclidean lattice Laplacian Q ■ The 
numerical value is v (0) « 0.16 0. 

Next, we consider the higher-order correction factor 
_^stress high-temperature expansion of the partition 

function Z — Zt^ooZ^^''^^^ ■ The correction factor is ob- 
tained by carrying out the integration over the displace- 
ment field u(x) in Eq. (|17|l . which makes (Ti2(x) discrete- 
valued, and further summing over the defect gauge fields 



^stress 



1^(1 -5x,e') 


E 




.cri2(x)e\/3Z/2_ 



+(5: 



x,B' 



X exp 



E 



E 



4/3A(l + i^) 



4a^ 

X i;a(x-x') — V(fe)e(fe)a;/V(„)e(„)j^/ 



0'12(x') 



(39) 



Most of the (712 configurations give a zero contribution to 
^■stress ^Yie sum Z](Ti2(x)eV3Z/2 EHl- We can extract 
the non-zero contributions by going to the gauge field 
x(x) with the definition cri2(x) = a^Va;Vyx(x). In order 
to get the summation values of the fields x(x) we have to 
determine the lowest value z such that za^V^^V^Jx.xo G 
\/?tL/2 for an arbitrary lattice position Xq. This value is 
given by z = 3/2. Neglecting boundary terms, we obtain 



^stress 



X exp 



n 



E 

x(x)e3Z 
1 



(40) 



4/3a (1 + 1^) 



5]x(x)(^.A)-^(x-x')x(x') 



where (x) is short for a"' (x) [recall 132(1 ] . Numer- 
ical integration yields 

(va)"'(0)«18.5 , (fA)"'(ae(,))«-4.43, (41) 

so that the lowest correction to the free energy is 

_^strcss ij^^j^ ^ 2Ne^^^^^^^^^^^^^^^^^^^'^^ . (42) 

The next term in the expansion of —(3F^'^°^^ is about a 

factor exp[— 9(ua)"^(0)/(16/3a(1 + '^))] smaller than the 
lowest-order term in 142() . By taking into account the 
lowest-order result /3a (1 + t^) ~ 0.6 we obtain from this 
an extremely small factor e^^". Thus, the lowest stress 
configuration of the high-temperature expansion in l(42() 
can be trusted for /3a(1 + J^) < 10. 

In the case of the square lattice one gets for Z^'^'^^^ 
the same formula as in 1(42(1 . except that the range of 
summation for the gauge field x(x) is X]x(x)gZ' ^^'^ that 
(ua)^^ has to be replaced by (wn)~^. The numerical 
values are (wn)~^(0) = 20 and (fn)~^(aej,,y) — —8. 

Using the above-derived results we can now calcu- 
late the melting temperature from the intersection of 
the low-temperature expansion Zt^o^a^ ^^'^ high- 
temperature expansion Zt^ooZ^'''°^^ of the partition 
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function Z. Instead of H28|l . we find from the corrections 
(|36|1 and H42I) that the melting point satisfies 



triangular lattice 



square lattice 



+ v) K, ^exp 
X exp 



. iVe '^^<i+'') 



(43) 



with /? — Pa in the case of the triangular lattice and 
f3 = /3n for the square lattice. The parameters are 





triangular 


square 


A 


0.51 


0.81 


Ml 


4 


4 


A'h 


12 




N 


4 


4 


Bi 


6.45 


6.31 


B2 


8.6 




C 


10.4 


5 



(44) 



Expression F/NksT = F/NksT + ln((l - + 
v))/N is a function of (3{\ + v). We show in Fig. 4 the 
corresponding curves of the low- and high-temperature 
expansion for the triangular lattice and the square lat- 
tice. We obtain for the square lattice that the high and 
low-temperature curves intersect in one point. In con- 
trast to this, the curves of the triangular lattice do not 
intersect. This is caused by an enhancement of the defect 
contributions to the low-temperature expansions for the 
triangular lattice. The two curves have merely a would- 
be intersection near 



/3a (1 



0.6, 



(45) 



where the distance between both curves on the /3(1 -I- v)- 
axis is approximatively 



A/3a(1 



0.05. 



(46) 



We point out that the non-intersecting of the two curves 
is not reasoned in a failure of our approximations. As 
estimated above the two curves can be trusted for 0.15 < 
I3a{1 + v) < 10. 

For a square lattice there was a definite intersection at 



/?□(!■ 



0.8. 



(47) 



In either case, the higher-order corrections to the melt- 
ing temperature are small compared to the lowest-order 
result H28|) . What is the reasons for the difference in the 
transition properties of the two lattices? We see from 
(|44|l that the numbers Mi and Bi arc almost identical for 
the triangular and the square lattice. The yy-configuration 
for the triangular lattice of this term is a triangle with 
charges ±2/3 shown in Fig. 3. The corresponding -q- 
configuration for the square lattice is given by a basic 
square with charges ±1 The difference between trian- 
gular and the square lattices lies in the term proportional 




0.2 0.4 0.6 0.8 I 1.2 1.4 1.6 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 

/3(l + i/) Pil + iy) 



FIG. 4: Low- and high-temperature expansions of F/NksT = 
F/NksT + ln((l - u^^^il + i^))/N for the triangular lattice 
and the square lattice. The thin line curves are the lowest- 
order results of Section V.A, while the thick line curves take 
into account also the defect and stress contributions of Sec- 
tion V.B. Note that stress corrections are negligible. The 
low-temperature lowest-order and higher-order curves for the 
triangular lattice are practically on top of each other. 



to A/2 in 143|) whose 77-configuration is shown in Fig. 2(a). 
A term of this type is absent for square lattices. The de- 
fects lie on two neighboring rhombuses of charges ±2/3. 
The corresponding configuration for the square lattice 
would consist of two neighboring squares. The distance 
between the rhombuses is much smaller than that of the 
two nearby squares. This results in an energy reduction 
of the corresponding defect configuration rj for the trian- 
gular lattice. This was the reason why this contribution 
was negligible for square lattices. 

From previous work we know that melting of two di- 
mensional crystals can proceed in only two possible sce- 
narios: either by the proliferation of defects as in three 
dimensions, in which case the transition is of first or- 
der, or by a two-step unbinding of defects (KTHNY). It 
was shown by Saito in Monte Carlo simulations of 
a gas of dislocations on a triangular lattice and in the 
textbook by general arguments that the core ener- 
gies of the dislocations decide which of the two options 
is chosen. For small core energies, the transition is of 
first order, for large core energies, the KTHNY scenario 
prevails. This agrees with the properties of our har- 
monic defect model. From the long-range behavior of 
the Coulomb force between the dislocations the unbind- 
ing transition will always take place even if the defects 
are not yet prolific at a temperature given by the univer- 
sal equation KT criterion. This implies that the /3-value 
at the melting point must always satisfy the inequality 
/3a (1 + v) > v^/tt « 0.55 irrespective of the lattice sym- 



metry . If the inequality is almost an equality, which 
is the case for our result (|45|l , the transition is of the 
KTHNY type. For square lattices, the situation is differ- 
ent. By comparing the melting temperature of the square 
(|47|1 and triangular H45() lattice we see that defects prolif- 
erate before dislocations can separate. This agrees with 
the naive expectation that square lattices should be less 
stable against thermal fluctuations than triangular ones. 
In the lattice models this is a consequence of the smaller 
core energies of dislocations on square lattices 0, as we 
know from Saito's simulation work. 

The behavior of the free energy shown in Fig. 4 ex- 
cludes the single-first-oder scenario for triangular lat- 
tices, leaving only the KTHNY scenario. Indeed, for this 
scenario we cannot expect an intersection of the curves in 
Fig. 4. The intermediate hexatic phase between the two 
continuous transitions, which the model should exhibit, 
would be inaccessible to both low- and high-temperature 
expansions. With less certainty, we extract an estimate 
for the separation of the two KT transitions from the dis- 
tance A/3 of the two nonintersecting curves. Our finding 
of a KTHNY scenario is in accordance with recent com- 
puter simulations for tri ang ular Lennard- Jones as well as 
electron lattices HHH HI ■ 

The situation for the square lattice is much clearer. 
There is a definite intersection point in Fig. 4, indicating 
that the melting transition is of first order, as found by 
computer simulations 0, . When trying to compare 
this result with real crystals in nature one has the dif- 
ficulty that square lattices are hard to produce. They 
need complicated interparticle forces, and only Weber et 
al. [iOI succeeded to do build such crystals with the help 
of a three-body potential in a computer model. They 
found indeed a first-order melting transition in their sim- 
ulation, in agreement with the prediction in 0, , and 
with the above conclusion. 



VI. LENNARD-JONES AND ELECTRON 
LATTICES IN 2D 

Let us compare our results quantitatively with com- 
puter simulations of 2D Lennard- Jones as well as electron 
lattices. In both cases, the ground state is triangular, 
due to the simplicity of the interaction potential. In the 
Lennard-Jones lattice, the interaction potential is 



14 J (r) = 4e 



(- 

V r 



12 



(48) 



Frenkel and McTargue |33 carried out an isothermal- 
isochoric molecular dynamics simulation, and observed 
a hexatic phase in accordance with KTHNY theory. In 
contrast, Abraham at. al. and others 33, 34] found 
with the help of both molecular dynamics as well as 
Monte-Carlo simulations that the melting transition of 
the Lennard-Jones lattice is of first order. The discrep- 
ancy induced simulations on larger systems, which found 
again KTHNY-like melting transitions [2lll2^ . 



Abraham 
Udink 

our calc. (w/o ren.) 
our calc. (w/ ren.) 




FIG. 5: Reduced melting temperature ksTm/t of Eq. 15111 as 
a function of the reduced density pa^. The dotted line denotes 
the melting temperature calculated with the T = elastic 
constants. The solid curve takes the thermal softening of 
the elastic constants by non-harmonic elasticity into account, 
calculated in Section VII. For comparison we show the mel ting 
temperature obtained by Abraham and Barker et al. 
and Udink et al. pH by simulation. The two curves of the 
simulation enclose the coexisting crystal-liquid region or the 
hexatic phase, respectively. 



To apply our model, we extract from the potential 148|) 
the elastic moduli A and /i as follows 2]: 



A = - 

V 



A* = - 

V 



180(-V^-36f-'^ 
a/ \a 



(49) 
(50) 



Inserting these into our result H45|l. we obtain for the 
melting temperature Tm of the triangular Lennard-Jones 
lattice: 



112(aV)^[(^V)'-0.3] 



1 



0.6 (27r) 



(51) 



4.725 (cr>)-^[(CT>)^-0.3]. 



Here p is the particle density 1/v. Figure 5 shows fcsTm/e 
as a function of cr^p in comparison with the two melting 
curves enclosing the coexisting crystal-liquid region or 
the hexatic phase, respectively, calculated by Abraham 
and Barker et al. in 1331 getting a first-order transition 
and by Udink et al. in |2l| consistent with KTHNY. The 
curves are calculated by determining the phase transi- 
tion line via thermodynamic integration. The results in 
|2lj | are in good accordance with the simulations on large 
systems showing also a KTHNY transition. 

Obviously, the value 0.6 in relation H45|) is too small 
to agree with the simulations. This discrepancy will be 
removed below in Section VII. 

Next, we discuss the system of a 2D electron lattice. 
Wigner [s^ predicted in 1934 that a gas of electrons 
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should becomes a triangular crystal at low temperature 
and density. Experimentally, this was first observed by 
Grimes and Adams 4] for electrons trapped on a liquid 
helium surface. The electrons interact via the Coulomb 
potential 



(52) 



and their kinetic energy is negligible at low temperature 
and density. As in the case of the Lennard-Jones lattice, 
the ground state of the system is built of a triangular 
lattice • Hockney and Brown !36j were the first which 
carry out a molecular dynamic simulation of this system. 
They found a A-like melt ing transition. In contrast, Gann 
et al. and Kalia et al. [33 found, via Monte-Carlo and 
molecular dynamics studies, a first-order melting tran- 
sition. Differing from both results, Morf 38] obtained 
a melting temperature consistent with the KTHNY the- 
ory, when he took into account the higher-order nonlinear 
corrections to the shear modulus. His results were re- 
produced bymore recent simulations on larger electron 
lattices plT2^. 

The elastic constants for the electron lattice are cal- 
culated in Refs. j2i 12] . The Lame constant is infinite, 
so that the electron lattice has an infinite modulus of 
compression. The shear constant is equal to 



^ = ^7eV^'^/^ ?? ~ 0.245. 



(53) 



The simulation data yield for the melting temperature 
Tm the relation 



r„ w 130 ± 10, 



(54) 



where F is the reduced inverse temperature 

r = ^e^/^keT. (55) 

To compare this result with our theory we insert the elas- 
tic constants into Eq. H45|l and find r„i ~ 86.4, which is 
about 30% smaller than the simulation value. Also here, 
we shall remove the discrepancy in Section VII. 

Concluding this section we observe that the interval 
A/3a(1 + v) k, 0.05 found in Fig. 4 for the triangular 
lattice in Eq. (|46|l agrees with the temperature interval 
of the intermediate (hexatic) phase seen in the computer 
simulations of both Lennard-Jones as well as electron lat- 
tices ill m mil. 



VII. ANHARMONIC CORRECTIONS TO 
MELTING FORMULA 

The melting temperature found in the last section from 
formula H45(l lies too high. The same problem was en- 
countered by Thouless [39j in his calculation of the melt- 
ing temperature for the electron lattice, who used the 
melting formula /3a (1 + i') — VS/^r ~ 0.55 where dis- 
locations unbind in the KTHNY theory. His formula is 
quite close to ours in Eq. 145|) . It was shown by Morf 



in [Sfll that the discrepancy between the prediction with 
the simulated melting temperature of the electron lattice 
can be explained by the thermal softening of the elas- 
tic constant /i, which he determined with the help of a 
computer simulation to be [recall (I5E 



Ai(r)«M0)(l-30.8/F), 



(56) 



where ^(0) is the zero-temperature shear modulus H53|l . 
By using H56|l one obtains from our formula (|45|l r,„ w 
117.2, which is now close to Morf's simulation value H54|l 
and in quite excellent agreement with the most recent 
melting value in Ref. 0| obtained by Monte-Carlo simu- 
lation of a large Wigner lattice F,„ « 123. 

For 2D Lennard-Jones crystals, the thermal softening 
can only be estimated theoretically. We obtain a rea- 
sonable approximation of the renormalized elastic con- 
stants by observing that the thermal softening due to 
non linear elasticity is much stronger for the transversal 
than for the longitudinal sound velocity 0. A further 
simplification comes from the fact 2] that the renormal- 
ization of the transversal sound is proportional to the 
square root of /x, and has in first approximation a uni- 
versal low-temperature law depending only on the lattice 
structure. This yields the approximations 



2^l{T) + A(T) w 2^(0) + A(0) , 
^(T) w /i(0) 



1 



/3(0) Vm(0) 2M0) + A(0) 



(57) 



(58) 



The constant c is the same for all lattices with the same 
structure. Its value can therefore be deduced from (|56|l : 



CA w 0.11 . 
Combining ()57|l and H58() we obtain 



u{T) 



1 



1.(0) -f (1 - Ko)) 



/3(0) 



HO) 



(59) 

(60) 
(61) 



Taking this renormalization into account, wc obtain from 
(I45II the melting formula for triangular lattices 



f3^{T)[l + y{T)) « /3A(0)(l + K0))-cK0)(3-i^(0)) 
« 0.6. (62) 

For the Lennard-Jones lattice where v{0) k, 0.48, this 
leads to a melting formula 



/3A(0)(l + i^(0))«0.74, 



(63) 



rather than H45() . thus lowering the melting tempera- 
ture by about 20%. The dashed curve in Fig. 5 shows 
the resulting melting temperature for the Lennard-Jones 
lattice, which is compatible with the simulation data 
'2l','33'|. It is obtained from (EH) by multiplying the right- 
hand side with the softening factor 0.6/0.74 w 0.81. 
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VIII. TWO-DIMENSIONAL LINDEMANN 
PARAMETER OF MELTING 

The Lindemann parameter expressed in terms 
of the expectation value of the lattice displacement 
\J {u^{Yi)) /a is not applicable in two dimensions, where 
the expectation value diverges due to excessive long- 
wavelength fluctuations of the displacement field. An 
ap pro priate modified Lindemann parameter in 2D is 



L = y ([-"(xj) - u(xj+i)] )/a , 



(64) 



where x^ and x^^.! are nearest neighbors in the lattice. 
The correlation function of the displacement field follow- 
ing from the energy (jHJ is 



K(k)uj(-k)> -7V(fcBT) (65) 




vfJ. g (i4r(,)if(„j)e(,)e(„))2 

On square lattices, we must merely replace 
(2/3)-f^(/)(e(;))a:,y ^ K^^y and (2/S)K^i){et^i))j;,y K^,y 
From this equation we obtain for both square and 
triangular lattices the modified Lindemann parameter 



(66) 



Let us compare our Lindemann parameter with that 
determined by Bedanov and Gadiyak via computer sim- 
ulations for Lennard- Jones and electron lattices ^3 ■ For 
the electron lattice we insert the temperature-dependent 
elastic constants fi{T) of Eq. (jSOJ and i/(r) = iy{0) ^ 1 
into the melting formula (|62|l . and obtain from H66|) 
La ~ 0.20. The corresponding simulation value of Be- 
danov and Gadiyak is La ~ 0.17, quite close to this. 

Now we turn to the Lennard- Jones lattice. Using the 
renormalized elastic constant (|61|l with 1/(0) « 0.48 and 
the melting formula 1)62(1 we obtain from ((66|l the Lin- 
demann parameter La ~ 0.20, again very close to the 
value La « 0.18 found in the simulations of Bedanov 
and Gadiyak [ioj . 

Hence our theory yields generalized Lindemann pa- 
rameters for both electron and Lennard- Jones lattices in 
good agreement with simulation results of Bedanov and 
Gadiyak 



a formula of the Lindemann type, with a modification due 
to the two dimensions. The value of the Lindemann pa- 
rameter is predicted, and agrees with estimates obtained 
from computer simulations on both Lennard- Jones and 
electron lattices. 

The melting transition is determined from the in- 
tersection of the free energies calculated once from 
a low-temperature expansion and once from a high- 
temperature expansion. While the square lattice melts 
in a first-order transition, the curves for the triangular 
lattice are compatible with two continuous transitions, as 
predicted by the KTHNY theory and found in computer 
simulations. The difference is due to the enhancement of 
the defect fluctuations in comparison to the square lat- 
tice caused by the smaller distance between basic nearby 
rhombuses on triangular lattices, in comparison to the 
corresponding defect configuration of two nearby squares 
on the square lattice. The resulting enhancement of 
smallest defect configurations leads to a decrease of the 
melting temperature. 



APPENDIX A: ELIMINATION OF THE GAUGE 
DEGREES OF FREEDOM 

In this section we eliminate the gauge degrees of free- 
dom in the n(^ijn) sum of Eq. H12|) enforced by the func- 
tional We shall prove the following: The gauge 
degrees of freedom in (|12|l with (|11|) and (|13|l are fixed in 
the case that one chooses the gauge fixed integer-valued 
defect fields (|15|l with the boundary condition H16|) for 
the bulk defect field n. 

By using and 113|) one can transform Hll|l to the 
following expression 

V:ru(x) -V(2)u(x) +n(2„)(x)e(,„) 



= -V(2)U- 



2^(21) - 71(22) 
^"(21) 



(Al) 



Vj;u(x) -> -^(2V(i) +V(2))u(x) 



V3 



(27i(im)e(„i) +n(2m)e(„)) 



1 



^(2V(i) +V(2))u(x) 



(A2) 



1 

71 



1(11) + 2n(i2) - 5'T-(21) + ^i(22) 

-V3n(ii) - 



Hii) ~ —"(21) 



IX. SUMMARY 

We have set up the simplest possible lattice model of 
defect melting on a two-dimensional triangular lattice. It 
accounts for the correct elastic fluctuations and by means 
of discrete- valued defect gauge fields for the fiuctuations 
of dislocations and disclinations. The latter give rise to 
melting transitions, and the melting temperature follows 



Here we put 7!(/3) = which is possible due to the over- 
counting of the basis e(/) . Now, we can carry out the sub- 
stitution of the displacement fields u(x) u(x)-t-aN(x) 
in (|A1|) where aN(x) is some displacement field corre- 
sponding to a jump from one lattice site to another. We 
determine it by the requirement that 

aV(2)iV:c(x) = n(2m)e(„)^ = ^'^(21) " "(22) (A3) 
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and 



^(2V(i)+V(2))A^,(x) 



- ■^(2"-(lm)e(„i)^ + n(2m)e(m)j^) - "(11) 



From these equations we obtain 



2^^(21) 



(A4) 



aV(,)7V,(x) e I , aV(,)iVy(x) e (A5) 



and 



(2a V(2)A^2:(x)) mod 2 



-^V(2)Af,(x) ) mod 2. 



(A6) 

To derive the last equation we have used the fact that 
the half-integer terms on the right hand sides of (|A3|) 
and ljA4|l coincide. We get unique solutions of ljA3|) and 
(jA4|) when fixing Nx{:>c) and Ny{x) on one half of the 
boundary of the system where we now suppose that we 
have approximately a square sample. Thus, we can for 
example fix the values of Nx{x) and Ny{x) on the up- 
per and the rightmost boundary of the sample which we 
denote by B. Then, due to the periodic boundary con- 
ditions for the displacement field u the values of N.j;{x.) 
and Ny (x) are determined on the whole boundary of the 
system. Because of the periodic boundary conditions for 
Nx{'x.) and Ny(x.) it is not clear that the unique solution 
of (|A3|) and (|A4|) respects these periodic boundary condi- 
tions. That this is not generally true can be most easily 
seen for the square lattice 2] . Whether this is true or not 
depends on the values of the defect configuration on the 



right hand sides of (|A3|) . ljA4p and the boundary values 
of 7V^(x) and Ny{x) on B. 

The values of the field n(x) in ()15|> or H17|l are now 
given by the following expression 



V3 
2 



a(V^iVj;(x) -f VyiV^(x)) 



+n(2m)e{m)y - "^(2"(lm) + "-(2m))e(m)^ 

= -2n(ii) + n(i2) 

+aV(i)iV^(x)-f aV3V(i)iVj,(x). (A7) 

Here we take into account that the Hamiltonian (TJ de- 
pends on the lattice derivates of the displacement field 
only in the strain combination V xUy -\- V yUx- We now 
observe that V(i) AT,, j,(x) can be written as a linear com- 
bination of V(2)iVa;^2^(x') with y' = y or y' = y + e(^i)^ and 
a lattice derivate corresponding to a rightmost boundary 
edge. From this, l|A5|l . (|A6|I . we obtain that n{x) on the 
whole lattice, given by the right-hand side of ljA7p . is 
integer valued when n(x) on the boundary x e S have 
integer values. 

The boundary values N^-ix) and Ny{'K) on B are fixed 
by choosing the boundary conditions H16() for n(x). The 
periodic boundary conditions for Nxi^x.) and Ny{'x.) on 
B can be fixed when taking into account IjAljl and (|A2|I 
only for s. ^ B. From lAip and l|A2p we obtain further 
that we can fix ri(22) (x) = Q on B. 

As a result we obtain the gauge-fixed path integral 114|l 
with (|15|l and (|16|l for the triangular lattice which takes 
into account all defect degrees of freedom. 
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